Counterion density profiles at charged flexible membranes 
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Counterion distributions at charged soft membranes are studied using perturbative analytical and 
simulation methods in both wealc coupling (mean-field or Poisson-Boltzmann) and strong coupling 
limits. The softer the membrane, the more smeared out the counterion density profile becomes and 
counterions pentrate through the mean-membrane surface location, in agreement with anomalous 
scattering results. Membrane-charge repulsion leads to a short-scale roughening of the membrane. 
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The study of charged colloids and biopolymers faces a 
fundamental problem: In theoretical investigations, the 
central object which is primarily computed is the charge 
density distribution in the electrolyte solution adjacent 
to the charged body 0. Experimentally measurable ob- 
servables are typically derived from this charge distribu- 
tion. For example, the force between charged particles 
follows from the ion density at the particle surfaces via 
the contact-value theorem. Likewise, the surface tension 
and surface potential are obtained as weighted integrals 
over the ion distributions. It has proven difficult to mea- 
sure the counterion distribution at a charged surface di- 
rectly because of the small scattering intensity. Notable 
exceptions are neutron scattering contrast variation with 
deuterated and protonated organic counterions |3| and 
local fluorescence studies on Zinc-ion distributions us- 
ing X-ray standing waves Clearly, direct comparison 
between theoretical and experimental ion distributions 
(rather than derived quantities) is desirable as it provides 
important hints how to improve theoretical modeling. 

In a landmark paper the problem of low scattering in- 
tensity was overcome by anomalous X-Ray scattering on 
stacks of highly charged bilayer membranes 01 • Anoma- 
lous scattering techniques allow to sensitively discrimi- 
nate counterion scattering from the background, and a 
multilayer consisting of thousands of charged layers gives 
rise to substantial scattering intensity. Since then, simi- 
lar techniques have been applied to charged biopolymers 
0, Q and to oriented charged bilayer stacks, where the 
problem of powder-averaging is avoided . 

However, scattering on soft bio-materials brings in a 
new complication, not considered theoretically so far: 
soft membranes and biopolymers fluctuate in shape, and 
thus perturb the counterion density proflle. Comparison 
with standard theories for rigid charged objects of sim- 
ple geometric shape becomes impossible. Here we fill this 
gap by considering the counterion-density profile close to 
a planar charged membrane which exhibits shape fluctu- 
ations governed by bending rigidity. As main result, we 
derive for a relatively stiff membrane closed-form expres- 
sions for the counterion density profile in the asymptotic 
low and high-charge limits which compare favorably with 



our simulation results. These parametric proflles, which 
exhibit a crucial dependence on the membrane stiffness, 
will facilitate the analysis of scattering results since they 
allow for a data fit with only a very few physical parame- 
ters. In previous experiments, a puzzling ion penetration 
into the lipid region was detected but interpreted as an 
artifact 4] . We show that ion penetration indeed occurs 
and is due to the correlated ion-membrane spatial fluc- 
tuations. The electrostatic coupling between membrane 
charges and counterions not only modifles the counte- 
rion density profile but also renormalizes the membrane 
roughness: the short-scale bending rigidity is reduced, 
charged membranes become locally softer. 

The Hamiltonian H — -I- i/e of the membrane- 
connterion system consists of the elastic membrane part 
Hm and the electrostatic part H,,. We discretize the 
membrane shape on a two-dimensional Nj^ x iV/, square 
lattice with lattice constant a and rescale all lengths by 
the Gouy-Chapman length /i — \l2'KqiBOm according to 
r = [IV, where Om = QMjN^c? is the projected charge 
density of the membrane and Ib = /AireoskBT is the 
Bjerrum length (e is the elementary charge, e the dielec- 
tric constant). Parametrizing the membrane shape by 
the height function /i(x), the elastic membrane energy in 
harmonic approximation reads in units of fc^T j^: 
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(1) 



where A is the Laplace operator, Kg is the bare bend- 
ing rigidity and g = gfi^ is the rescaled strength of the 
harmonic potential. The electrostatic energy accounts 
for the interaction of N counter-ions of valence q and 
M membrane charges of valence Q, related by the elec- 
troneutrality condition QM = qN, 
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where s — 2'!Tq igam denotes the coupling parame- 
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FIG. 1: Simulation snapshots for a) H = 0.2, ^^/m = 0.80, 
Ko = 0.07, g = 0.57, a = 0.18, d = 2.2 and b) H = 1000, 
Cx/m = 0.38, ifo = 174, 3 = 0.006, 5 = 13.21, d = 160. 
The simulations are done using TV — 100 counter-ions and 
M — 100 membrane-ions on a Nl ~ 60 x 60 membrane lattice. 

ter 1^. The rescaled position of the ith counterion is 
fj while the A:-th membrane-ion is located at R^. = 
(xfc,/i(xfe) — d) where the membrane charges are dis- 
placed hy d = 2aNLM~^/'^ beneath the membrane sur- 
face which is impenetrable to the point-like counterions. 
This way we can largely neglect charge-discreteness ef- 
fects [lO| and concentrate on shape-fluctuation effects. 
In most of our simulations the membrane ions are mobile 
and move freely on the membrane lattice, with a packing 
fraction C = M/iV|. For the long-ranged electrostatic 
interactions we employ laterally periodic boundary con- 
ditions using Lekner-Sperb methods P|. To minimize 
discretization and finite-size effects, the number of lat- 
tice sites Nl and the rescaled strength of the harmonic 
potential g are chosen such that the lateral height-height 
correlation length of the membrane obeys the inequal- 
ity: a < = (4A'o/5)^/^ < iVi,a|3. Simulations are run 
for typically 10^ Monte Carlo steps using 100 counter- 
ions and 100 membrane ions. In Fig^we show two sim- 
ulation snapshots. The counter-ions form in the weak 
coupling hmit (S — 0.2, Fig^Ja) a diffuse dense cloud 
while in the strong coupling limit (S = 1000, Fig^Jb, 
note the anisotropic rescaling) the lateral ion-ion dis- 
tances are large compared to the mean separation from 
the membrane. Pronounced correlations between mem- 
brane shape fluctuations and counterion positions are ob- 
served in both snapshots. 

The qualitatively different ionic structures at low/high 
coupling strength are reflected by fundamentally different 



analytic approaches in these two limits: Starting point is 
the exact expression for the partition function 




By performing a Hubbard-Stratonovich transformation 
and a transformation to the grand-canonical ensemble, 
we arrive at the partition function [lH : 

Z ~ J VhVcjje-""^^''^-"'^^''-'^-''^/^ (4) 

The field i(j> is the fluctuating electrostatic potential [lll |. 
The electrostatic action reads 

/ - Hi)) (5) 

where 9{z) — 1 for z > and zero otherwise. The ex- 
pectation value of the counter-ion density is calculated 
by the help of the generating field 7r(r) according to 
(p(f)) ~ 2^TE5\nZ/S^^{r)fJ,^^ and reads 

(^(f)) = :#%-^(^(^-Mx))e-^«) . (6) 

The dimensionless fugacity A is determined by the 
normalization condition of the counterion distribution 
/ dr (p(r)) — N, which is in rescaled units equivalent to 

A / df (^e{S - /i(i))e-'"*(*')^ = 1. The partition function 
Eq.0 is intractable. In the weak coupling limit, S — > 0, 
fluctuations of the field (f) around the saddle point value 
are small and gaussian variational methods become ac- 
curate The variational Gibbs free energy reads: 

= Fo + (H^[h, ct>, 7r]/S + H„,[h] - Ha[h, q}]) ^ (7) 

Here (• • • )o is an average with the variational hamiltonian 
Hq and _Fo is the corresponding free energy. The most 
general Gaussian variational hamiltonian Hq is 

Ho[h,(j)] = ^ J dxdx/h{x)K^^{x,x/)h{x/) 

+ ^ J drdr/n{r)v-^{r,r/)n{r/) , (8) 

where the field is defined by fl{r) := 0(f) — (j)o{r) + 
i J dx/dx//P(f;x/)ii'"i(x/,x//)/i(x//) and P is the con- 
nected correlation function P(f; x/) = (i(/)(f )/i(x/))g. The 
variational parameters are the mean potential 0o, the 
coupling operator P, the propagator of the electrostatic 
field V and the membrane propagator K. For K we 
use the bare propagator of the uncharged membrane 
K{ii,x/) = -4(|J[)2kei(\/2|x - x/|/|[?)/7r, where the 
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bare membrane roughness ^5. is given by l/y/GAKog = 
(Ci.)^ — (^^(0))o 1^1' Assuming the charge propagator v 
to be isotropic and translational invariant (which is an 
approximation) v turns out to be the bare Coulomb prop- 
agator, v{r) = 1/r. The remaining variational equations 
6Fy/6P — 5Fy/5(t)Q — are solved perturbatively in an 
asymptotic small ^5. expansion, i.e. for a relatively stiff 
membrane. The solution for P for x = x/ is expressed 
in terms of the Meijer's Q function and reads (neglecting 
terms of 0{{i1_f )): 
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The result for the mean potential (f>Q is given by 
Ea. (|10|l and reduces in the limit ^ to the known 
Gouy-Chapmann potential i(j>{z) = 21n(l + zlfisl 
We defined the auxiliary function w{z) as: w{z) :~ 



J 



2(^0)2/^exp{~fV2(Cl)2}-5erfc(z7^2(e0)2). The 
counterion density is calculated according to Eq.® and 
up to third order in ^± given by Eq. Hll() ; it reduces to the 
known mean-field counter-ion density {p(z)) = (1 -I- zY^ 
in the case of vanishing membrane roughness ■ 
In FigI21we show the laterally averaged counterion den- 
sity profiles for weak coupling S — 0.2 obtained from MC 
simulation (solid squares) for several membrane rough- 
nesses . For the comparison with the analytical expres- 
sion Ea. Hll() (solid lines) we use the discrete membrane 
propagator K:^^ — AKf){cofi\^T:n / N l] + cos[27rm/iVi] — 
2)"^ /a^A-g and calculate the membrane roughness accord- 
ing to {£!\_Y = n ^rnn- The lateral correlation length 
follows as — l/(2^'[5^/^). For z > ^° the counterion 
profile approaches the corresponding profile for a planar 
surface, but for z < we find pronounced deviations 
from the flat surface profile. For — 1.211 the analyti- 
cal result and the simulation result disagree, showing the 
limitation of our small ^5. expansion. 



i(j)Q{z) 



w{~z) + 21n [l + z - ~zw{~z)/A - (IV2)'erf (7^)] + 0{{£,lf) : 5 > 



25 - ^2 + yj(^z) (1 - z/2) - [S^lfeii 



/2 + o{{e^y)-rz<o 



(10) 



-i(t>o(z) 




1 - ^erf 



r 



-2P±(5)- 



/2^^« 



\3\ 



(11) 



In the strong coupling limit S — > cx) we expand the parti- 
tion function Q in inverse powers of S j9| . Starting point 
is the exact expression Eq.JBJ. After some manipulation 
we find for the leading term: 
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^ f drl5[zl-h{xl))v(r,rl) _|_(9(2~1). (12) 



This strong coupling expansion is equivalent to a 
virial expansion, and hence the leading term cor- 
responds to the interaction of a single counterion 
with a fiuctuating charged membrane 9]. For stiff 
membranes we can employ a small-gradient expansion, 
^ / df/(5(z/-/i(x/)i>(f-f/) ~ C'--z-|-/(if//i(x/)/^(f,f/), 
where C is an unimportant constant and the func- 
tion /^(f) is defined by: /^(f,f/) := S{zf — h{xf))x 

V 3/2 



(^(5 - z/)-(i-i/) • V/h{i/)j/2n (|i-i/|2 + (5-z/)2)' 
Expanding Ea. l|12|l in powers of gives rise to: 
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The density H13|) reduces to the known SC density 



(p(z)) = in the limit ^5. ^0 We compare in 
Fig|21the analytically obtained counterion density profiles 
(solid lines) with the laterally averaged densities obtained 
using MC simulations (open triangles) for S = 1000 and 
different • The analytic approximation reproduces the 
simulated profiles very well. Similar to the weak cou- 
pling case, the profiles approach the corresponding strong 
coupling density for counter-ions at a planar surface for 
z 3> ^5., but deviate noticeable from the planar distri- 
bution for z < Comparison of mobile and immo- 
bile membrane ions gives no detectable difference for the 
counterion profic (Fig. 2 inset). 

In the analytics so far we used the bare membrane 
roughness without modification due to electrostatics. 
In FigOwe show the ratio of the membrane rough- 
ness measured in the MC simulation, and for the 
bare uncharged membrane, as a function of the coupling 
parameter S for two different surface fractions C (open 
symbols). The ratio is larger than unity, i.e. charges on 
the membrane increase the roughness. This short-range 
roughening, which allows membrane charges to increase 
their mutual distance and is thus not area-preserving, 
has to be distinguished from the electrostatic stiffening 
in the long-wavelength limit which has been predicted 
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FIG. 2: Rescaled counterion density {p{z)) = {p{z)) /27r^_ga^ 
as a function of the rescaled distance z = z/jj. from Monte 
Carlo simulations (data points) and asymptotic theory (solid 
lines). In the weak coupling limit (H = 0.2, solid squares), 
the membrane roughness is = 1.211,0.3184,0 and ^jj — 
0.2483 , 0.2933 , oo from bottom to top. In the strong cou- 
pling limit (H = 1000, open triangles) we have ~ 
1.211,0.3184,0 and = 17.2475 , 20.7458 , oo from bottom 
to top. Numerical errors are smaller then the symbol sizes. 
In all cases the membrane-ions are mobile and the packing 
fraction is ^ = 0.028. The inset compares profiles for H — 0.2, 
i1 = 0.3184 for C = 0.028 (diamonds) and C = 0.25 (circles) 
for mobile membrane ions and resuhs for H = 0.2, = 1.211, 
— 0.028 for mobile (squares) and fixed (stars) membrane 
ions and H = 1000, il = 1.211, ( = 0.028 for mobile (trian- 
gle) and fixed (crosses) membrane ions. 

on the mean- field-level El E [13. Local roughening 
corresponds to protrusion degrees of freedom of single 
lipids. Yet a distinct softening mechanism, effective at 
intermediate wave leng ths, is due to electrostatic corre- 
lations effects 0, Il9l Eoj , which is missed by standard 
mean-field approaches. Experimentally, both membrane 
stiffening |2h | and, for highly charged membranes, soften- 
ing has been observed 22]. To distinguish effects due to 
membrane charges and counterions we calculate via ex- 
act enumeration and within harmonic approximation the 
membrane propagator Kmn for a charged discrete mem- 
brane without counterions. The roughness ratio from 
this analytical calculations is shown as a solid line, and 
again cross-checked by MC simulations without counte- 
rions (filled symbols) . The good agreement with the MC 
data containing counterions shows that the softening ef- 
fect is mostly due to the repulsion of charges on the mem- 



brane itself. Experimentally, this short-scale roughening 
will show up in diffuse X-ray scattering data. 
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FIG. 3: Ratio of simulated and bare roughness £,±/£,± as a 
function of H for — 0.028 and — 0.3184 (open squares) 
and 1° = 1.2111 (open stars), C = 0.25 and I" = 0.3184 
(open triangles) and — 1.2111 (open diamonds). The solid 
lines and solid symbols are analytical and MC results without 
counterions ((" = 0.028 lower branch, ^ = 0.25 upper branch). 
The inset shows the ratio as a function of the packing 

fraction ( for H = 0.2 (squares) and H = 1000 (triangles), 
= 0.3184 in both cases. 



